home *** CD-ROM | disk | FTP | other *** search
/ SGI Hot Mix 17 / Hot Mix 17.iso / HM17_SGI / research / lib / ts_smooth.pro < prev    next >
Text File  |  1997-07-08  |  6KB  |  145 lines
  1. ; $Id: ts_smooth.pro,v 1.5 1997/01/15 03:11:50 ali Exp $
  2. ;
  3. ; Copyright (c) 1995-1997, Research Systems, Inc.  All rights reserved.
  4. ;       Unauthorized reproduction prohibited.
  5. ;+
  6. ; NAME:
  7. ;       TS_SMOOTH
  8. ;
  9. ; PURPOSE:
  10. ;       This function computes central, backward or forward moving-averages
  11. ;       of an n-element time-series (X). Autoregressive forecasting and 
  12. ;       backcasting is used to extrapolate the time-series and compute a 
  13. ;       moving-average for each point of the time-series. The result is an 
  14. ;       n-element vector whose type is identical to X.
  15. ;
  16. ; CATEGORY:
  17. ;       Statistics.
  18. ;
  19. ; CALLING SEQUENCE:
  20. ;       Result = TS_SMOOTH(X, Nvalues)
  21. ;
  22. ; INPUTS:
  23. ;       X:    An n-element vector of type float or double containing time-
  24. ;             series samples, where n >= 11.
  25. ;
  26. ; Nvalues:    A scalar of type integer or long integer that specifies the
  27. ;             number of time-series values used to compute each moving-average.
  28. ;             If central-moving-averages are computed (the default), this
  29. ;             parameter must be an odd integer of 3 or greater.
  30. ;
  31. ; KEYWORD PARAMETERS:
  32. ;   BACKWARD: If set to a non-zero value, backward-moving-averages are 
  33. ;             computed. The Nvalues parameter must be an integer greater 
  34. ;             than 1.
  35. ;
  36. ;     DOUBLE: If set to a non-zero value, computations are done in
  37. ;             double precision arithmetic.
  38. ;
  39. ;    FORWARD: If set to a non-zero value, forward-moving-averages are computed.
  40. ;             The Nvalues parameter must be an integer greater than 1.
  41. ;
  42. ;      ORDER: A scalar of type integer or long integer that specifies 
  43. ;             the order of the autoregressive model used to compute the 
  44. ;             forecasts and backcasts of the time-series. Central-moving-
  45. ;             averages require Nvalues/2 forecasts and Nvalues/2 backcasts.
  46. ;             Backward-moving-averages require Nvalues-1 backcasts.
  47. ;             Forward-moving-averages require Nvalues-1 forecasts. 
  48. ;             A time-series with a length in the interval [11, 219] will use
  49. ;             an autoregressive model with an order of 10. A time-series with
  50. ;             a length greater than 219 will use an autoregressive model with 
  51. ;             an order equal to 5% of its length. The ORDER keyword is used to
  52. ;             override this default.
  53. ;
  54. ; EXAMPLE:
  55. ;       Define an n-element vector of time-series samples.
  56. ;         x = [6.63, 6.59, 6.46, 6.49, 6.45, 6.41, 6.38, 6.26, 6.09, 5.99, $
  57. ;              5.92, 5.93, 5.83, 5.82, 5.95, 5.91, 5.81, 5.64, 5.51, 5.31, $
  58. ;              5.36, 5.17, 5.07, 4.97, 5.00, 5.01, 4.85, 4.79, 4.73, 4.76]
  59. ;
  60. ;       Compute the 11-point central-moving-averages of the time-series.
  61. ;         result = ts_smooth(x, 11)
  62. ;
  63. ; REFERENCE:
  64. ;       The Analysis of Time Series, An Introduction (Fourth Edition)
  65. ;       Chapman and Hall
  66. ;       ISBN 0-412-31820-2
  67. ;
  68. ; MODIFICATION HISTORY:
  69. ;       Written by:  GGS, RSI, September 1995
  70. ;       Modified:    GGS, RSI, July 1996
  71. ;                    Modified keyword checking and use of Double precision.
  72. ;-
  73.  
  74. FUNCTION TS_Smooth, x, Nvalues, Backward = Backward, Double = Double, $
  75.                                 Forward = Forward, NaN = NaN, Order = Order
  76.  
  77.   ON_ERROR, 2
  78.  
  79.   TypeX = SIZE(x)
  80.   nX = TypeX[TypeX[0]+2]
  81.  
  82.   ;Check time-series length.
  83.   if nX lt 11 then $
  84.     MESSAGE, "Time-series input must be a vector of at least 11 elements."
  85.  
  86.   ;If the DOUBLE keyword is not set then the internal precision and
  87.   ;result are identical to the type of input.
  88.   if N_ELEMENTS(Double) eq 0 then Double = (TypeX[TypeX[0]+1] eq 5)
  89.  
  90.   ;Define output type.
  91.   if Double eq 0 then SMx = FLTARR(nX) else SMx = DBLARR(nX)
  92.  
  93.   ;Order of Autoregressive model.
  94.   ;A time-series with a length in the interval [11, 219] will use an 
  95.   ;Autoregressive model with an Order of 10. A time-series with a length
  96.   ;greater than 219 will use an Autoregressive model with an Order equal 
  97.   ;to 5% of its length. The ORDER keyword is used to override this default.
  98.   if KEYWORD_SET(Order) eq 0 then $
  99.     Order = MAX([10L, LONG(0.05 * nX)])
  100.  
  101.   if KEYWORD_SET(Backward) ne 0 then begin 
  102.     ;Compute Backward-moving-averages.
  103.     if Nvalues lt 2 then MESSAGE, $ ;Nvalues must be 2 or greater.
  104.       "Backward average; Nvalues must be an integer greater than 1."
  105.     ;Requires (Nvalues-1) backcasted values.
  106.     x = [TS_FCAST(x, Order, Nvalues-1, /BACKCAST, Double = Double), x]
  107.     for k = 0, nX-1 do $
  108.       SMx[k] = TOTAL(x[k:Nvalues+(k-1)], Double = Double)
  109.     ;Restore x to its input state.
  110.     if TypeX[TypeX[0]+1] eq 4 then x = FLOAT(x[Nvalues-1:*]) $
  111.     else x = DOUBLE(x[Nvalues-1:*])
  112.     RETURN, SMx/Nvalues
  113.   endif else if KEYWORD_SET(Forward) ne 0 then begin 
  114.     ;Compute Forward-moving-averages. 
  115.     if Nvalues lt 2 then MESSAGE, $ ;Nvalues must be 2 or greater.
  116.       "Forward average; Nvalues must be an integer greater than 1."
  117.     ;Requires (Nvalues-1) forecasted values.
  118.     x = [x, TS_FCAST(x, Order, Nvalues-1, Double = Double)]
  119.     for k = 0, nX-1 do $
  120.       SMx[k] = TOTAL(x[k:Nvalues+(k-1)], Double = Double)
  121.     ;Restore x to its input state.
  122.     if TypeX[TypeX[0]+1] eq 4 then x = FLOAT(x[0:nX-1]) $
  123.     else x = DOUBLE(x[0:nX-1])
  124.     RETURN, SMx/Nvalues
  125.   endif else begin
  126.     ;Compute central-moving-averages.
  127.     if Nvalues lt 3 then MESSAGE, $ ;Nvalues must be odd and 3 or greater.
  128.       "Central average; Nvalues must be an odd integer greater than 2."
  129.     if LONG(Nvalues) MOD 2 eq 0 then SMwidth = LONG(Nvalues) + 1 $
  130.     else SMwidth = LONG(Nvalues)
  131.     ;Requires (Nvalues/2) forecasted values and (Nvalues/2) backcasted values;
  132.     ;where Nvalues is an odd integer.
  133.     x = [TS_FCAST(x, Order, SMwidth/2, /BACKCAST, Double = Double), x, $
  134.          TS_FCAST(x, Order, SMwidth/2, Double = Double)]
  135.     for k = 0, nX-1 do $
  136.       SMx[k] = TOTAL(x[k:SMwidth+(k-1)], Double = Double)
  137.     ;Restore x to its input state.
  138.     if TypeX[TypeX[0]+1] eq 4 then x = FLOAT(x[SMwidth/2:SMwidth/2+(nX-1)]) $
  139.     else x = DOUBLE(x[SMwidth/2:SMwidth/2+(nX-1)])
  140.     RETURN, SMx/SMwidth    
  141.     ;SMx = smooth(ts, SMwidth) & SMx = SMx[SMwidth/2:SMwidth/2+nX-1]
  142.   endelse
  143.  
  144. END
  145.